function result = DynesFit_swave_gauss_fraction_Base (Param, eV) %Param(1) = delta0 (meV), Param(2)=gamma, Param(3)=T (mK), Param(4) = slopeBG, Param(5) = offsetBG, Param(6) = sigma

delta0 = Param(1);
sigma = Param(5);
w = 1; % fraction of SC region
Temperature = 116.4; %Temperature of the sample (mK)
beta = 1/(8.61733E-5 .* Temperature);  % 8.61733E-5 is kb/e
xrange = -10/beta : 0.2/beta: +10/beta; % Energy integration section for Fermi-Dirac distribution

if delta0-4*sigma < 0
    yrange = 0 : 4*sigma/50 : delta0+4*sigma; % Delta integration section for gaussian distribution at E>0
    yrange2 = delta0-4*sigma: 4*sigma/50 : 0; % Delta integration section for gaussian distribution at E<0
    [x,delta] = meshgrid(xrange,yrange);
    [x2,delta2] = meshgrid(xrange,yrange2);

    Gauss = 1./(sqrt(2*pi).*sigma).*exp(-(delta-delta0).^2/(2*sigma.^2)); %Gaussian distribution for E>0
    Gauss2 = 1./(sqrt(2*pi).*sigma).*exp(-(delta2-delta0).^2/(2*sigma.^2)); %Gaussian distribution for E<0

    %Dynes = @(E) abs(real( (E-1i.*Param(2)) ./ sqrt((E-1i.*Param(2)).^2-delta.^2)))*Gauss(delta);  % Here energy also has the unit meV
    % FD = beta .* exp(beta.*(xrange))./ ((exp(beta.*(xrange))+1.0).^2); % Fermi-Dirac distribution

    innerFtn1 = (1-w) .* beta .* exp(beta.*(xrange))./ ((exp(beta.*(xrange))+1.0).^2); % Normal state contribution
    innerFtn2 = w .* Gauss2.* beta .* exp(beta.*(x2))./ ((exp(beta.*(x2))+1.0).^2); % State not SC yet contribution (E<0 in Gaussian)
    innerFtn3 = @(E) w .* abs(real( (x+E-1i.*Param(2)) ./ sqrt((x+E-1i.*Param(2)).^2-delta.^2))).*Gauss.* beta .* exp(beta.*(x))./ ((exp(beta.*(x))+1.0).^2); % State already SC contribution (E>0 in Gaussian)

    result =  cell2mat(arrayfun(@(E) ((Param(3).*E+ Param(4)).*( trapz(xrange,innerFtn1) + trapz(xrange,trapz(yrange2,innerFtn2,1)) + trapz(xrange,trapz(yrange,innerFtn3(E),1)))),eV,'UniformOutput',false));

elseif delta0-4*sigma > 0
    yrange = delta0-4*sigma : 4*sigma/50 : delta0+4*sigma; % Delta integration section for gaussian distribution at E>0
    [x,delta] = meshgrid(xrange,yrange);

    Gauss = 1./(sqrt(2*pi).*sigma).*exp(-(delta-delta0).^2/(2*sigma.^2)); %Gaussian distribution for E>0

    innerFtn1 = (1-w) .* beta .* exp(beta.*(xrange))./ ((exp(beta.*(xrange))+1.0).^2); % Normal state contribution
    innerFtn3 = @(E) w .* abs(real( (x+E-1i.*Param(2)) ./ sqrt((x+E-1i.*Param(2)).^2-delta.^2))).*Gauss.* beta .* exp(beta.*(x))./ ((exp(beta.*(x))+1.0).^2); % State already SC contribution (E>0 in Gaussian)

    result =  cell2mat(arrayfun(@(E) ((Param(3).*E+ Param(4)).*( trapz(xrange,innerFtn1) + trapz(xrange,trapz(yrange,innerFtn3(E),1)))),eV,'UniformOutput',false));
end

end

